library(RUVSeq); library(EDASeq); library(RColorBrewer)
colData(se)$mix=NULL; colData(se)$sf=NULL
ruv <- RUVg(round(zfGenes), spikes, k=1)
differences <- matrix(data=c(c(1,2,4,5), c(3,6,7,-1)), byrow=TRUE, nrow=2)
ruvs <- RUVs(round(zfGenes), genes, k=1, differences)
require(edgeR)
design <- model.matrix(~strain, data=colData(se))
y <- DGEList(counts=assays(se)$counts[genes,], group=colData(se)$strain)
y <- calcNormFactors(y, method="upperquartile")
y <- estimateGLMCommonDisp(y, design)
y <- estimateGLMTagwiseDisp(y, design)
fit <- glmFit(y, design)
res <- residuals(fit, type="deviance")
ruvr <- RUVr(round(zfGenes), genes, k=1, res)
myHeatmap <- function(d, ...) {
results(d) %>% .[order(.$padj),] %>% head(50) %>% rownames %>%
counts(d)[.,] %>% '+'(1) %>% log2 %>%
annHeatmap(ann=colData(d), cluster=list(cuth=15, col=colors), ...) %>% plot
}
colors <- brewer.pal(3, "Set2")